library(aqp)
library(soilDB)
library(farver)
library(cluster)
library(ape)
library(colorspace)
library(latticeExtra)
library(knitr)
library(MASS)
library(vegan)
Compute several indices of color contrast from 2 vectors of Munsell colors, m1 and m2. Comparisons fully vectorized, e.g. comparisons performed element-wise \(m1_i ~\leftrightarrow~ m2_i\) \(...\) \(m1_n ~\leftrightarrow~ m2_n\).
Results as a data.frame.
# examples
m1 <- c('10YR 6/3', '7.5YR 3/3', '10YR 2/2', '7.5YR 3/4')
m2 <- c('5YR 3/4', '7.5YR 4/4', '2.5YR 2/2', '7.5YR 6/3')
# result is a data.frame
d <- colorContrast(m1, m2)
# check
kable(d, row.names=FALSE)
| m1 | m2 | dH | dV | dC | dE00 | cc |
|---|---|---|---|---|---|---|
| 10YR 6/3 | 5YR 3/4 | 2 | 3 | 1 | 31.112562 | Prominent |
| 7.5YR 3/3 | 7.5YR 4/4 | 0 | 1 | 1 | 9.522489 | Faint |
| 10YR 2/2 | 2.5YR 2/2 | 3 | 0 | 0 | 6.959582 | Faint |
| 7.5YR 3/4 | 7.5YR 6/3 | 0 | 3 | 1 | 30.054654 | Distinct |
Results summarized in a figure.
# more examples
m1 <- c('10YR 6/3', '7.5YR 3/3', '10YR 2/2', '7.5YR 3/4', '2.5Y 6/8', '5B 4/6')
m2 <- c('5YR 3/4', '7.5YR 4/4', '2.5YR 2/2', '7.5YR 6/3', '10YR 2/1', '5GY 3/4')
# graphical comparison
colorContrastPlot(m1, m2)
Simulated redoximorphic feature colors, contrast classes and \(\Delta{E_{00}}\).
m1 <- paste0('7.5YR 4/', 2:8)
m2 <- rep('10YR 5/2', times=length(m1))
colorContrastPlot(m1, m2, labels = c('F3M', 'MAT'), d.cex = 0.8, col.cex = 0.8)
m1 <- paste0('5Y 4/', 5:1)
m2 <- rep('10YR 3/5', times=length(m1))
colorContrastPlot(m1, m2, labels = c('RMX/FED', 'MAT'), d.cex = 0.8, col.cex = 0.8)
m1 <- c('2.5Y 3/2', '5Y 3/2', '2.5GY 3/2', '5GY 3/2', '7.5BG 3/2', '2.5B 3/2')
m2 <- rep('10YR 3/5', times=length(m1))
colorContrastPlot(m1, m2, labels = c('RMX', 'MAT'), d.cex = 0.8, col.cex = 0.8)
Test rules and cases outlined in Soil Survey Technical Note 2.
10YR 6/3 vs 5YR 3/4
contrastClass(v1=6, c1=3, v2=3, c2=4, dH=2, dV=3, dC=1, verbose = TRUE)
## $faint
## v1 c1 v2 c2 dH dV dC f.case1 f.case2 f.case3 low.value.chroma res
## 1 6 3 3 4 2 3 1 FALSE FALSE FALSE FALSE Prominent
##
## $distinct
## v1 c1 v2 c2 dH dV dC d.case1 d.case2 d.case3 res
## 1 6 3 3 4 2 3 1 FALSE FALSE FALSE Prominent
7.5YR 3/3 vs 7.5YR 4/4
contrastClass(v1=3, c1=3, v2=4, c2=4, dH=0, dV=1, dC=1, verbose = TRUE)
## $faint
## v1 c1 v2 c2 dH dV dC f.case1 f.case2 f.case3 low.value.chroma res
## 1 3 3 4 4 0 1 1 TRUE FALSE FALSE FALSE Faint
##
## $distinct
## v1 c1 v2 c2 dH dV dC d.case1 d.case2 d.case3 res
## 1 3 3 4 4 0 1 1 FALSE FALSE FALSE Faint
10YR 2/2 vs 2.5YR 2/2
contrastClass(v1=2, c1=2, v2=2, c2=2, dH=0, dV=0, dC=0, verbose = TRUE)
## $faint
## v1 c1 v2 c2 dH dV dC f.case1 f.case2 f.case3 low.value.chroma res
## 1 2 2 2 2 0 0 0 TRUE FALSE FALSE TRUE Faint
##
## $distinct
## v1 c1 v2 c2 dH dV dC d.case1 d.case2 d.case3 res
## 1 2 2 2 2 0 0 0 FALSE FALSE FALSE Faint
7.5YR 3/4 vs 7.5YR 6/3
contrastClass(v1=3, c1=4, v2=5, c2=3, dH=0, dV=3, dC=1, verbose = TRUE)
## $faint
## v1 c1 v2 c2 dH dV dC f.case1 f.case2 f.case3 low.value.chroma res
## 1 3 4 5 3 0 3 1 FALSE FALSE FALSE FALSE Distinct
##
## $distinct
## v1 c1 v2 c2 dH dV dC d.case1 d.case2 d.case3 res
## 1 3 4 5 3 0 3 1 TRUE FALSE FALSE Distinct
Differences in Munsell hue (\(dH\)) are computed according to radial position (clock-wise) from 5R   to    5PB.
Use huePosition() to extract hues in order, or convert a vector of hues into an index of positions.
# hues used in the description of soil color (except for neutral hues)
hues <- huePosition(x=NULL, returnHues = TRUE)
# convert hue names into position
hue.pos <- huePosition(hues)
# check
kable(head(cbind(hue.pos, hues)))
| hue.pos | hues |
|---|---|
| 1 | 5R |
| 2 | 7.5R |
| 3 | 10R |
| 4 | 2.5YR |
| 5 | 5YR |
| 6 | 7.5YR |
Arranged in CIELAB colorspace, the ordering of hues looks like this. The huePositionPlot function in {sharpshootR} will make this kind of figure for any combination of value and chroma.
# make these into real colors by fixing value and chroma at '6'
colors <- paste0(hues, ' 6/6')
# convenient labels for figure
hue.labels <- sprintf("%s\n%s", hues, hue.pos)
# combine colors names, hex representation of colors, and CIELAB coordinates
x <- data.frame(
colors,
hex=parseMunsell(colors),
parseMunsell(colors, returnLAB=TRUE),
stringsAsFactors = FALSE
)
# make the figure
plot(B ~ A, data=x, type='n', asp=0.75, xlab='A-coordinate', ylab='B-coordinate', ylim=c(-25, 45), main='Hue Order per TN #2\nCIELAB Colorspace')
abline(h=0, v=0, lty=2)
points(B ~ A, data=x, pch=15, col=x$hex, cex=4.5)
text(x$A, x$B, labels = hue.labels, cex=0.75)
Inspect the data used to generate the figures with the additional argument returnData=TRUE.
Basic chart, single reference hue.
contrastChart(m = '7.5YR 4/3', hues = c('7.5YR'))
Threshold color chips by \(\Delta{E_{00}}\).
contrastChart(m = '7.5YR 4/3', hues = c('7.5YR'), thresh = 15)
Multiple reference hues.
contrastChart(m = '7.5YR 4/3', hues = c('5YR', '7.5YR', '10YR'))
Split chips by contrast class, single reference hue.
contrastChart(m = '7.5YR 4/3', hues = c('7.5YR'), style = 'CC')
Split chips by contrast class, multiple reference hues. This works better when run through latticeExtra::useOuterStrips().
fig <- contrastChart(m = '7.5YR 4/3', hues = c('5YR', '7.5YR'), style = 'CC')
useOuterStrips(fig, strip = strip.custom(bg=grey(0.85)), strip.left = strip.custom(bg=grey(0.85)) )
# load full LUT from {aqp}
data(munsell)
# subset some common colors
hues <- c('10R', '2.5YR', '5YR', '7.5YR', '10YR', '2.5Y', '5Y')
x <- subset(munsell, subset=value %in% 3:5 & chroma %in% 1:3 & hue %in% hues)
# convert into an ordered factor, according to hue position
huePos <- huePosition(x=NULL, returnHues = TRUE)
x$hue <- factor(x$hue, levels=huePos, ordered = TRUE)
x$hue <- droplevels(x$hue)
# ok
table(x$hue)
##
## 10R 2.5YR 5YR 7.5YR 10YR 2.5Y 5Y
## 9 9 9 9 9 9 9
nrow(x)
## [1] 63
# convert into hex notation for plotting
x$color <- munsell2rgb(x$hue, x$value, x$chroma)
## Notice: converting hue to character
x$munsell <- sprintf("%s %s/%s", x$hue, x$value, x$chroma)
Arrange colors in a familiar manner.
xyplot(value ~ chroma | hue, data=x,
xlim=c(0.5, 3.5), ylim=c(2.5, 5.5),
scales=list(alternating=3, x=list(at=1:3), y=list(at=3:5)),
as.table=TRUE, strip=strip.custom(bg='grey'),
subscripts=TRUE,
panel=function(xx, yy, subscripts, ...) {
panel.grid(-1, -1)
panel.points(xx, yy, col=x$color[subscripts], pch=15, cex=6)
}
)
Evaluate \(\Delta~E_{00}\) via {farver} pacakge.
# dE00
# pair-wise distances are only present in the upper triangle
# lower triangle is 0
d <- compare_colour(x[, c('L', 'A', 'B')], from_space='lab', white_from = 'D65', method='cie2000')
# copy color codes for convenience
dimnames(d) <- list(x$munsell, x$munsell)
# convert to dist object, which expects distances in the lower triangle
d.dist <- as.dist(t(d))
# hierarchical clustering
h <- as.phylo(as.hclust(diana(d.dist)))
# nMDS
mds <- MASS::sammon(d.dist, trace = FALSE)
# remove 0's in lower triangle and diagonal
d[lower.tri(d, diag = TRUE)] <- NA
hist(d)
# double-check that ordering is correct
# yes
table(h$tip.label == x$munsell)
##
## TRUE
## 63
par(mar=c(0.5,0.5,1,0.5), bg='black', fg='white')
plot(h, label.offset=0.55, edge.color='white', tip.color='white', cex=0.66, direction='downwards', font=1)
tiplabels(pch=15, cex=1.5, col=x$color, offset=0.23)
title('Divisive Hierarchical Clustering', col.main='white', line=0)
mtext(text = expression("Distance metric:"~Delta*E['00']), side = 3, adj = 0, line=-0.5)
plot(h, edge.color='white', tip.color='white', cex=0.5, direction='downwards', show.tip.label=FALSE)
tiplabels(pch=15, cex=1.5, col=x$color, offset=0.25)
tiplabels(text=x$value, cex=0.55, frame='none', col='white', offset=0.23)
title('Divisive Hierarchical Clustering', col.main='white', line=0)
mtext(text = expression("Distance metric:"~Delta*E['00']), side = 3, adj = 0, line=-0.5)
par(mar=c(0.5,0.5,1.5,0.5), bg='black', fg='white')
plot(mds$points, type='n', axes=FALSE)
abline(h=0, v=0, col=grey(0.85), lty=3)
points(mds$points, pch=15, col=x$color, cex=3.5)
text(mds$points[, 1], mds$points[, 2], x$value, cex=0.75)
mtext(text = expression("Distance metric:"~Delta*E['00']), side = 3, adj = 0)
mtext(text = 'Munsell value', side = 1, adj = 1, line=-1)
title('nMDS', col.main='white', line=0.75)
plot(mds$points, type='n', axes=FALSE)
abline(h=0, v=0, col=grey(0.85), lty=3)
points(mds$points, pch=15, col=x$color, cex=3.5)
text(mds$points[, 1], mds$points[, 2], paste0(x$value, '/', x$chroma), cex=0.66)
mtext(text = expression(Delta*E['00']), side = 3, adj = 0)
mtext(text = 'Munsell value/chroma', side = 1, adj = 1, line=-1)
title('nMDS', col.main='white', line=0.75)
plot(mds$points, type='n', axes=FALSE)
abline(h=0, v=0, col=grey(0.85), lty=3)
points(mds$points, pch=15, col=x$color, cex=3.5)
text(mds$points[, 1], mds$points[, 2], x$hue, cex=0.66)
mtext(text = expression(Delta*E['00']), side = 3, adj = 0)
mtext(text = 'Munsell hue', side = 1, adj = 1, line=-1)
title('nMDS', col.main='white', line=0.75)
plot(mds$points, type='n', axes=FALSE)
abline(h=0, v=0, col=grey(0.85), lty=3)
points(mds$points, pch=15, col=x$color, cex=3.5)
text(mds$points[, 1], mds$points[, 2], x$munsell, cex=0.66)
mtext(text = expression(Delta*E['00']), side = 3, adj = 0)
mtext(text = 'Munsell color', side = 1, adj = 1, line=-1)
title('nMDS', col.main='white', line=0.75)
# locate color pairs which are at the edge of perceptable differences
idx <- which(d < 4)
row.idx <- row(d)[idx]
col.idx <- col(d)[idx]
low.dE00 <- cbind(dimnames(d)[[1]][row.idx], dimnames(d)[[2]][col.idx])
# not so helpful viz of pairs
swatchplot(list(
A=parseMunsell(low.dE00[, 1]),
B=parseMunsell(low.dE00[, 2])
))
# compute median delta-E00 by Munsell chip
l <- list()
for(i in seq_along(x$munsell)) {
m <- x$munsell[i]
d.i <- na.omit(c(d[i, ], d[, i]))
l[[i]] <- data.frame(m=m, stat=median(d.i), stringsAsFactors = FALSE)
}
# flatten
E00.summary <- do.call('rbind', l)
# check
kable(E00.summary[grep('10YR', E00.summary$m), ], digits=2)
| m | stat | |
|---|---|---|
| 10 | 10YR 3/1 | 11.15 |
| 11 | 10YR 3/2 | 10.64 |
| 12 | 10YR 3/3 | 11.92 |
| 13 | 10YR 4/1 | 10.45 |
| 14 | 10YR 4/2 | 10.35 |
| 15 | 10YR 4/3 | 11.37 |
| 16 | 10YR 5/1 | 11.82 |
| 17 | 10YR 5/2 | 11.58 |
| 18 | 10YR 5/3 | 13.32 |
# combine chip summary with source data
z <- merge(x, E00.summary, by.x='munsell', by.y='m', all.x=TRUE)
# chip size is proportional to median delta-E00
xyplot(value ~ chroma | hue, data=z,
xlim=c(0.5, 3.5), ylim=c(2.5, 5.5),
scales=list(alternating=3),
as.table=TRUE, strip=strip.custom(bg='grey'),
subscripts=TRUE,
panel=function(xx, yy, subscripts, ...) {
panel.grid(-1, -1)
panel.points(xx, yy, col=z$color[subscripts], pch=15, cex=z$stat * 0.4)
}
)
# subset some common colors
hues <- c('7.5YR')
x <- subset(munsell, subset=value %in% 3:8 & chroma %in% c(1,2,3,4,6,8) & hue %in% hues)
# convert into hex notation for plotting
x$color <- munsell2rgb(x$hue, x$value, x$chroma)
x$munsell <- sprintf("%s/%s", x$value, x$chroma)
# dE00
# pair-wise distances are only present in the upper triangle
# lower triangle is 0
d <- compare_colour(x[, c('L', 'A', 'B')], from_space='lab', white_from = 'D65', method='cie2000')
# copy color codes for convenience
dimnames(d) <- list(x$munsell, x$munsell)
# convert to dist object, which expects distances in the lower triangle
d.dist <- as.dist(t(d))
# nMDS
mds <- MASS::sammon(d.dist, trace = FALSE)
# rotate 90 degrees CCW
# to roughly follow Munsell color book page layout
# column-order
m <- matrix(
c(0, -1,
1, 0),
byrow = FALSE, ncol = 2
)
# apply transformation
mds.trans <- mds$points %*% m
par(mar=c(2.5, 0.5, 3, 0.5), bg='black', fg='white')
plot(mds.trans, type='n', axes=FALSE, asp=1)
abline(h=seq(-25, 25, 5), v=seq(-25, 25, 5), col='white', lty=3)
points(mds.trans, pch=22, bg=x$color, cex=6)
text(mds.trans[, 1], mds.trans[, 2], paste0(x$value, '/', x$chroma), cex=0.66)
# mtext(text = expression(Delta*E['00']), side = 3, adj = 0)
# mtext(text = '7.5YR', side = 1, adj = 1, line=-1)
title('Perceptual Distances\n7.5YR Page', col.main='white', line=0.75)
mtext(text = expression(Delta*E['00']%->%nMDS%->%90*degree~CCW~rotation), side = 1, adj = 0)
Map perceptual arrangement to original representation in color book.
p <- procrustes(mds.trans, x[, c('value', 'chroma')], scale = TRUE, symmetric = TRUE)
par(mar=c(2.5, 0, 3, 0), bg='black', fg='white')
plot(p, pch=22, kind=0, axes=FALSE, xlab='', ylab='')
points(p, display='rotated', pch=22, bg=x$color, cex=5)
points(p, display='target', pch=21, bg=x$color, cex=4)
lines(p, type='arrows', length=0.1)
text(p, display='rotated', cex=0.66, labels=paste0(x$value, '/', x$chroma))
title('Perceptual vs Color Book Arrangement', col.main='white', line=1)
title(sub='7.5YR Page', col.sub='white', line=0)
legend('bottomright', legend=c('color book', 'perceptual'), pch=c(22, 21), pt.bg=parseMunsell('7.5YR 4/6'), pt.cex=2, bty='n', inset=c(0.25, 0))
x <- fetchOSD('musick', colorState = 'dry')
y <- fetchOSD('musick', colorState = 'moist')
m1 <- sprintf("%s %s/%s", x$hue, x$value, x$chroma)
m2 <- sprintf("%s %s/%s", y$hue, y$value, y$chroma)
cc <- colorContrast(m1, m2)
colorContrastPlot(m1, m2, labels=c('Dry', 'Moist'), d.cex = 0.9)
Tinker with plot settings and optimize for dark background.
x <- fetchOSD('Leefield', colorState = 'dry')
y <- fetchOSD('Leefield', colorState = 'moist')
m1 <- sprintf("%s %s/%s", x$hue, x$value, x$chroma)
m2 <- sprintf("%s %s/%s", y$hue, y$value, y$chroma)
par(bg='black', fg='white')
colorContrastPlot(m1, m2, labels=c('Dry', 'Moist'), printMetrics = TRUE, d.cex = 1.25, col.cex = 1.5, label.cex = 2, label.font = 2)
# fresh start
data(munsell)
# for now, most commonly used hue pages, but expanded to full range in value / chroma
x <- subset(munsell, subset = hue %in% c('10R', '2.5YR', '5YR', '7.5YR', '10YR', '2.5Y', '5Y'))
# create all pair-wise combinations
cols <- paste0(x$hue, ' ', x$value, '/', x$chroma)
z <- combn(cols, 2)
# convert to 2-column data.frame
# ~ 674,541 combinations!
z <- data.frame(t(z), stringsAsFactors = FALSE)
## TODO: make this a new function
# borrowed from aqp::colorContrast
m1.pieces <- parseMunsell(z[[1]], convertColors = FALSE)
m2.pieces <- parseMunsell(z[[2]], convertColors = FALSE)
# convert to value and chroma to numeric
m1.pieces[[2]] <- as.numeric(m1.pieces[[2]])
m1.pieces[[3]] <- as.numeric(m1.pieces[[3]])
m2.pieces[[2]] <- as.numeric(m2.pieces[[2]])
m2.pieces[[3]] <- as.numeric(m2.pieces[[3]])
# difference in number of hue chips, clock-wise, as specified in:
# https://www.nrcs.usda.gov/wps/portal/nrcs/detail/soils/ref/?cid=nrcs142p2_053569
dH <- abs(huePosition(m1.pieces[[1]]) - huePosition(m2.pieces[[1]]))
# difference in number of value chips
dV <- abs(m1.pieces[[2]] - m2.pieces[[2]])
# difference in number of chroma chips
dC <- abs(m1.pieces[[3]] - m2.pieces[[3]])
## ^^^^ need a new function
# just adjacent chips
# ~ 11,498 combinations
idx <- which(dH %in% c(0,1) & dV %in% c(0,1) & dC %in% c(0,1))
z <- z[idx, ]
# compute color contrast metrics on adjacent chips
d <- colorContrast(z[[1]], z[[2]])
# summaries
quantile(d$dE00)
## 0% 25% 50% 75% 100%
## 0.000000 4.765694 7.994458 9.797715 17.252535
# encode combinations
d$code <- sprintf("%s-%s-%s", d$dH, d$dV, d$dC)
# compute select quantiles of dE00 by combination
adj <- lapply(split(d, d$code), function(i) {
qq <- round(quantile(i$dE00, probs = c(0.05, 0.5, 0.95)))
data.frame(
dH = i$dH[1],
dV = i$dV[1],
dC = i$dC[1],
q05 = qq[1],
q50 = qq[2],
q95 = qq[3]
)
}
)
adj <- do.call('rbind', adj)
# all combinations
knitr::kable(adj)
| dH | dV | dC | q05 | q50 | q95 | |
|---|---|---|---|---|---|---|
| 0-0-1 | 0 | 0 | 1 | 0 | 2 | 4 |
| 0-1-0 | 0 | 1 | 0 | 6 | 8 | 10 |
| 0-1-1 | 0 | 1 | 1 | 7 | 9 | 11 |
| 1-0-0 | 1 | 0 | 0 | 1 | 4 | 6 |
| 1-0-1 | 1 | 0 | 1 | 2 | 5 | 7 |
| 1-1-0 | 1 | 1 | 0 | 6 | 10 | 13 |
| 1-1-1 | 1 | 1 | 1 | 6 | 10 | 13 |
# just changes along a single dimension
knitr::kable(adj[row.names(adj) %in% c('0-0-1', '1-0-0', '0-1-0'), ])
| dH | dV | dC | q05 | q50 | q95 | |
|---|---|---|---|---|---|---|
| 0-0-1 | 0 | 0 | 1 | 0 | 2 | 4 |
| 0-1-0 | 0 | 1 | 0 | 6 | 8 | 10 |
| 1-0-0 | 1 | 0 | 0 | 1 | 4 | 6 |
library(rms)
library(tactile)
data(munsell)
x <- subset(munsell, subset=value %in% 3:8 & chroma %in% 1:6 & hue %in% c('10R', '2.5YR', '5YR', '7.5YR', '10YR', '2.5Y', '5Y'))
cols <- paste0(x$hue, ' ', x$value, '/', x$chroma)
z <- combn(cols, 2)
z <- data.frame(t(z), stringsAsFactors = FALSE)
d <- colorContrast(z[[1]], z[[2]])
kable(head(d), row.names = FALSE)
| m1 | m2 | dH | dV | dC | dE00 | cc |
|---|---|---|---|---|---|---|
| 10R 3/1 | 10R 3/2 | 0 | 0 | 1 | 5.266520 | Faint |
| 10R 3/1 | 10R 3/3 | 0 | 0 | 2 | 9.026377 | Distinct |
| 10R 3/1 | 10R 3/4 | 0 | 0 | 3 | 11.849272 | Distinct |
| 10R 3/1 | 10R 3/5 | 0 | 0 | 4 | 14.044517 | Prominent |
| 10R 3/1 | 10R 3/6 | 0 | 0 | 5 | 15.819853 | Prominent |
| 10R 3/1 | 10R 4/1 | 0 | 1 | 0 | 8.811914 | Faint |
bwplot(cc ~ dE00, data=d,
par.settings=tactile.theme(),
xlab=expression(Delta*E['00']),
main='CIE Delta-E00 vs. Soil Color Contrast Class\nhue 10R-2.5Y, value 3-8, chroma 1-6',
scales=list(x=list(tick.number=10)),
panel=function(...) {
panel.grid(h=3, v=-1)
panel.bwplot(...)
})
bwplot(dV ~ dE00, data=d,
subset = dH < 1 & dC < 1,
par.settings=tactile.theme(),
xlab=expression(Delta*E['00']),
main='CIE Delta-E00 vs. Change in Munsell Value\nhue 10R-2.5Y, value 3-8, chroma 1-6',
sub = 'Constant Hue and Chroma',
scales=list(x=list(tick.number=10)),
panel=function(...) {
panel.grid(h=3, v=-1)
panel.bwplot(...)
})
bwplot(dH ~ dE00, data=d,
subset = dV < 1 & dC < 1,
par.settings=tactile.theme(),
xlab=expression(Delta*E['00']),
main='CIE Delta-E00 vs. Change in Munsell Hue\nhue 10R-2.5Y, value 3-8, chroma 1-6',
sub = 'Constant Value and Chroma',
scales=list(x=list(tick.number=10)),
panel=function(...) {
panel.grid(h=3, v=-1)
panel.bwplot(...)
})
bwplot(dC ~ dE00, data=d,
subset = dV < 1 & dH < 1,
par.settings=tactile.theme(),
xlab=expression(Delta*E['00']),
main='CIE Delta-E00 vs. Change in Munsell Chroma\nhue 10R-2.5Y, value 3-8, chroma 1-6',
sub = 'Constant Hue and Value',
scales=list(x=list(tick.number=10)),
panel=function(...) {
panel.grid(h=3, v=-1)
panel.bwplot(...)
})
idx.F <- which(d$cc == 'Faint' & d$dE00 > 12)
idx.D <- which(d$cc == 'Distinct' & d$dE00 > 20)
idx.P <- which(d$cc == 'Prominent' & d$dE00 < 20)
kable(head(d[idx.F, ]), row.names = FALSE)
| m1 | m2 | dH | dV | dC | dE00 | cc |
|---|---|---|---|---|---|---|
| 10R 3/1 | 10R 5/1 | 0 | 2 | 0 | 18.68596 | Faint |
| 10R 3/1 | 10R 5/2 | 0 | 2 | 1 | 18.94273 | Faint |
| 10R 3/2 | 10R 5/1 | 0 | 2 | 1 | 19.83064 | Faint |
| 10R 3/2 | 10R 5/2 | 0 | 2 | 0 | 18.72647 | Faint |
| 10R 3/2 | 10R 5/3 | 0 | 2 | 1 | 18.78994 | Faint |
| 10R 3/2 | 5Y 3/2 | 6 | 0 | 0 | 13.85876 | Faint |
# plot(dE00 ~ dH, data=d)
# plot(dE00 ~ dV, data=d)
# plot(dE00 ~ dC, data=d)
dd <- datadist(d)
options(datadist="dd")
# there is a little bit of curvature
(m <- ols(dE00 ~ rcs(dH, 3) + rcs(dV, 3) + rcs(dC, 3), data=d))
## Linear Regression Model
##
## ols(formula = dE00 ~ rcs(dH, 3) + rcs(dV, 3) + rcs(dC, 3), data = d)
##
## Model Likelihood Discrimination
## Ratio Test Indexes
## Obs 31626 LR chi2 100532.26 R2 0.958
## sigma2.3709 d.f. 6 R2 adj 0.958
## d.f. 31619 Pr(> chi2) 0.0000 g 12.552
##
## Residuals
##
## Min 1Q Median 3Q Max
## -9.01452 -1.78228 -0.08724 1.68269 10.22195
##
##
## Coef S.E. t Pr(>|t|)
## Intercept 6.1933 0.0490 126.37 <0.0001
## dH 0.4303 0.0236 18.24 <0.0001
## dH' 1.1100 0.0325 34.14 <0.0001
## dV 4.4707 0.0260 171.90 <0.0001
## dV' 3.7974 0.0290 130.99 <0.0001
## dC 0.4921 0.0260 18.92 <0.0001
## dC' 0.7916 0.0290 27.31 <0.0001
##
# (m <- ols(dE00 ~ dH + dV + dC, data=d))
plot(Predict(m, dH=NA, dV=0, dC=0))
plot(Predict(m, dH=0, dV=NA, dC=0))
plot(Predict(m, dH=0, dV=0, dC=NA))
plot(Predict(m, dH=NA, dV=0:2, dC=0))
plot(Predict(m, dV=NA, dH=0:2, dC=0))
plot(Predict(m, dC=NA, dV=0:2, dH=0))
plot(summary(m, dH=c(0,1), dV=c(0,1), dC=c(0,1)))
plot(summary(m))
plot(anova(m), what='partial R2')
plot(nomogram(m))
For next time, how do you use these tables?
print(nomogram(m, dH=c(0, 1, 2, 3)))
## Points per unit of linear predictor: 2.535368
## Linear predictor units per point : 0.39442
##
##
## dH Points
## 0 0
## 1 1
## 2 3
## 3 6
##
##
## dV Points
## 0 0
## 1 12
## 2 27
## 3 49
## 4 74
## 5 100
##
##
## dC Points
## 0 0
## 1 1
## 2 3
## 3 7
## 4 11
## 5 15
This document is based on aqp version 1.29.